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the two- loop amplitude at the production threshold s > 4m 2 . 
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1. Introduction 

The top quark is the heaviest fermion of the Standard Model. By studying its properties 
in detail, it is hoped to elucidate the origin of particle masses and the mechanism of 
electroweak symmetry breaking. Since its discovery at the Fermilab Tevatron [1] a little 
more than a decade ago, its mass has been measured to within a few per cent, while its 
production cross section and couplings are only known with larger uncertainty. With the 
large number of top quarks expected to be produced at the LHC, the study of its properties 
will become precision physics. To interpret these upcoming precision data, equally precise 
theoretical predictions are mandatory. These demand foremost the calculation of higher 
order corrections in perturbative QCD. 

At present, the top quark pair production cross section [2-8] is known to next-to- 
leading order (NLO) in the QCD coupling constant, the same precision is available for 
single top production [9], tt+jet production [10] and tt+Z-boson production [11]. For the 
pair production cross section, resummation [12-14] of logarithmically enhanced corrections 
(next-to-leading logarithm, NLL) to all orders in the coupling constant improves upon the 
fixed-order NLO prediction. Electroweak one-loop corrections to ti production are equally 
available [15, 16]. 

Especially for the top quark pair production cross section, which is expected to be mea- 
sured to within a few per cent accuracy, it is believed that the current NLO+NLL prediction 
is not yet sufficiently accurate. Detailed recent studies [17-19] indicate a scale uncertainty 
on these predictions of 7%, and a parton distribution uncertainty of 6%. While the latter 
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Figure 1: Tree-level amplitude. Massive quarks are indicated by a thick line. 

may be improved upon by more precise determinations of the parton distribution functions 
in view of recent and upcoming data from HERA and LHC, the former requires the cal- 
culation of perturbative corrections at next-to-next-to-leading order (NNLO) in QCD. By 
approximating these corrections with the fixed-order expansion of the NLL prediction, one 
finds [17] a projected NNLO scale uncertainty of 3%, which is below the parton distribution 
uncertainty, and in line with the anticipated experimental error. 

The calculation of the full NNLO corrections to the top quark pair production cross 
section requires three types of ingredients: two- loop matrix elements for qq — > ti and 
gg —* ti, one-loop matrix elements for hadronic production of tt+(l parton) and tree-level 
matrix elements for hadronic production of ti+{2 partons). The latter two ingredients 
were computed previously in the context of the NLO corrections to ii+jet production [10]. 
They contribute to the ti production cross section through configurations where up to two 
final state partons can be unresolved (collinear or soft), and their implementation thus may 
require further developments of subtraction techniques at NNLO. 

Both two-loop matrix elements were computed analytically in the small-mass expansion 
limit s, \t\, \u\ ^> m 2 in [20,21], starting from the previously known massless two-loop matrix 
elements for qq — > q'q' [22] and gg — > qq [23]. An exact numerical representation of the 
two-loop matrix element qq — > ti has been obtained very recently [24]. It is the aim of the 
present paper to compute all two-loop contributions to qq — > ti arising from closed fermion 
loops in a compact analytic form, which provide a first independent validation of the recent 
results of [20,24], allow for a fast numerical evaluation, and permit the analytical study of 
the behavior of the top quark production cross section at threshold. 

This paper is structured as follows. In Section ^, we define our notation and kinematical 
conventions. Sections [|| and || describe the details of the calculation of the two-loop integrals 
and of the renormalization of the amplitudes. The results are presented and discussed 
in Section ||. We enclose two appendices describing the special functions used in our 
calculation and documenting the newly computed master integrals. 

2. Notation and Conventions 

We consider the scattering process 

q(pi) + qfa) — » tips) + *(P4) , (2.1) 
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in Euclidean kinematics, where pf = for i = 1, 2 and p| = —m 2 for i = 3,4. The 
Mandelstam variables are defined as follows 



Ol + P2) 2 , t = -(pi- p 3 ) 2 , u = - (pi - p 4 )' 



(2.2) 



Conservation of momentum implies that s + i + u = 2m 2 . 

The squared tree-level matrix element (averaged over the spin and color of the incoming 
quarks and summed over the spin of the outgoing ones), calculated in d = 4—2e dimensions, 
can be expanded in powers of the strong coupling constant as as follows: 



\M\ 2 (s,t,m,e) 



A7r 2 a 2 
N 2 



Ao + (— ) -4 , + 



IT 



(2.3) 



The tree-level amplitude involves a single diagram (Fig. |Tj) and its contribution to Eq. (^ 
is given by 



A = 4N C C F 



(t - m 2 ) 2 + {u- m 2 ) 2 2m 2 
1 — tt 1 — H e 



(2.4) 



where N c is the number of colors and Cf = (iV^ — 1)/2N C . As it is well known, the term 
proportional to the dimensional regulator e in Eq. ( |2.4| ) is mass independent. 

The NLO term A\ in Eq. ( |2.3D arises from the interference of one-loop diagrams with 
the tree-level amplitude [2-8]. The NNLO term A2 consists of two parts, the interference 
of two-loop diagrams with the Born amplitude and the interference of one-loop diagrams 
among themselves: 



The latter term ^l^ 1 ^ was studied extensively in [25]. A^ x0 \ originating from the two- 
loop diagrams, can be decomposed according to color and flavor structures as follows: 



(lxl) 



A 



(2x0) 



N C C F 



C 



N 2 



E, 



NfA + B + — + N t { N c Di + -±)+N h [ N c D t + 



El 

N c 



+N?Fi + NiN h F lh + NlF h 



(2.5) 



where Ni and Nh are the number of light- and heavy-quark flavors, respectively. The 
coefficients A, B, . . . , in Eq. ( |2.5[) are functions of s, t, and m, as well as of the dimen- 
sional regulator e. Recently, these quantities were calculated in [20] in the approximation 
s, \t\, \u\ 3> m 2 . For a fully differential description of top quark pair production at NNLO, 

(2x0) 

the complete mass dependence of A 2 is required. An exact numerical expression for 
it has been obtained very recently in [24]. In this work, we provide independent confir- 
mations of the recent results of [20, 24] by deriving exact analytic expressions for all the 
terms in Eq. ( p.5|) arising from two-loop diagrams involving at least a fermion loop (i.e. 
the coefficients Di,E{, Fj with i = l,h and j = l,h, Ih). 
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Figure 2: Non-reducible topologies for the light quark corrections. 



3. Calculation 



The calculation starts from the two-loop Feynman diagrams for qq —* tt, which are gen- 
erated using QGRAF [26], interfered with the tree-level amplitude, and simplified using 
FORM [27]. Out of the 218 two-loop diagrams contributing to the amplitude, 28 are pro- 
portional to Ni, 29 are proportional to N^, 2 are proportional to iVjiV^, while just one 
contributes to the Nf and N% parts. Most importantly, there is only one two- loop box 
topology contributing to the iVj part of the squared amplitude, and a single other two-loop 
box topology proportional to N^. These two box topologies are very similar to the ones en- 
countered in the evaluation of the two-loop QED corrections to Bhabha scattering [28,29], 
and can be evaluated with the same techniques. 

All two-loop integrals appearing in these amplitudes are reduced to a set of master 
integrals (Mis) by means of the standard method based on the Laporta algorithm [30]. Only 
part of these Mis were available in the literature [31-37] from previous two-loop calculations 
of the heavy quark form factors [38] and amplitudes for Bhabha scattering [28,29,39]. For 
the remaining integrals, we employed the differential equation method [40]. 

The reduction to Mis was carried out with two independent implementations of the 
Laporta algorithm, and large parts of it were cross checked with the Maple package A. I .R. 
[41]. The 12 irreducible topologies encountered in the calculation of the diagrams with 
a light quark loop are shown in Fig. |2[ The diagrams proportional to Nh also contain 
12 irreducible topologies, which can be found in Fig. ||. In both figures, thick internal 
lines indicate massive propagators, while thin lines indicate massless ones. An external 
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Figure 3: Non-reducible topologies for the heavy quark corrections. 



dashed leg carries a squared momentum s; other external lines indicate particles on their 
mass-shell, where pf = for thin lines and pf = —m 2 for thick lines. 

The analytic expressions of the one-loop Mis in Fig. |2| and Fig. || are well known. 
The large majority of the two-loop Mis is also known in the literature: explicit analytic 
expressions for all the two-loop Mis with the exception of the ones belonging to topologies 
Fig. |-(k), Fig. g-(l), Fig. g-(k), and Fig. can be found in [31-33,35]. 

The Mis associated with topologies Fig. g-(k), Fig. @-(l), Fig. g-(k), and Fig. §-(1) 
that were not available in the literature are collected in Appendix [b[ In calculating the 
Mis by means of the differential equation method, it is crucial to fix the undetermined 
integration constant(s) appearing in the solution of the differential equations. While there 
is no general method available to fix such initial condition, it is usually sufficient to know 
the behavior of the MI in some particular kinematic point; for example, knowing that the 
integral is regular for a certain value of s, one can impose the regularity of the solution of 
the differential equation in that point. This can be sufficient to determine the integration 
constant. In our calculation, the initial conditions for the single master integral belonging 
to topology Fig. §-(k) and the two Mis belonging to topology Fig. ||-(1) were determined by 
imposing the regularity of the solution of the differential equation in t = 0. However, this is 
not always sufficient. For topology Fig. |3|-(k), which has two Mis, imposing the regularity 
of both Mis in t = allowed to fix only one of the two initial conditions required. In 
order to fix the second integration constant, we had to use another piece of information, 
namely that the scalar integral with all the denominators raised to power one diverges 
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at most logarithmically at the threshold t = m 2 . The final result for these two Mis was 
then checked by calculating their t — > limit with the Mellin Barnes technique, using 
the Mathematica packages Ambre [42] and MB [43]. For what concerns the MI of the box 
topology in Fig. |2|-(1), the initial conditions were defined calculating the integral in t = 
with Mellin Barnes techniques. 

All the Mis were calculated in the non-physical region s < 0, where they are real and 
can be conveniently written as functions of the dimensionless variables 

V— s + 4m 2 — y/—s —t —u 

V-s + 4m 2 + V^s m m 

The Mis of the topology Fig. |3|-(e) are an exception. In this case it is convenient to employ 

the variable 

v^i-V- S -4m 2 (32) 



\J — s + v 7 — s — 4m 2 

The transcendental functions appearing in the Mis are one- and two-dimensional harmonic 
polylogarithms (HPLs) [44,45]. In the result one finds one-dimensional HPLs of maximum 
weight four and two-dimensional HPLs of maximum weight three. Both sets of functions 
can be rewritten in terms of conventional Nielsen's polylogarithms. In Appendix we 
briefly review the definition of the HPLs employed and we collect the expression of some 
of them in terms of Nielsen's polylogarithms. 

Following the procedure outlined in the present section, it is possible to obtain the 
expression of the bare squared matrix elements involving diagrams proportional to Ni 
and/or N^. After this goal is achieved, it is then necessary to renormalize the ultraviolet 
divergencies. In the next section, we briefly discuss the renormalization procedure and we 
explicitly list the needed renormalization constants. 

4. Renormalization 

The renormalized QCD matrix element is obtained from the bare one by expanding the 
following expression : 

A-en = JJ Zj£ >n .4-bare (<%,bare — * Z as as > "^bare — ¥ Z rn,mj , (4.1) 

n 

where Z WY ,n is the external leg wave function renormalization factor, as is the renormalized 
coupling constant and m is the renormalized heavy quark mass. (In the rest of the section 
we suppress the subscript "5"' in as)- 

We postpone the discussion of mass renormalization to the end of the section and we 
start by considering the coupling constant and wave function renormalization. 

We introduce the following quantities: 

7T 7T 

By expanding the amplitude and the wave function renormalization factor in do we find: 

Aen(«bare) = a oA) + OqAx + alA 2 + O(a ) , 

Z WF<n = 1 + a Q 6z£l n + ajjtfzS, n + O(og) . (4.3) 
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The relation between ao and a is given by: 

a = a + a 2 5Z^ + a 3^(2) + ( a 4) _ (44) 



By employing Eqs. (|4.3,4.4) in Eq. ( |4.l| ) we find 

Aen = aA + a 2 4en + a 3 A$L + 0{a A ) , 



4en = Ai + \*Z { An + SZ^j Ao , 
A^ n = A 2 +(^2 ^wl,n + 25ZW j ^ + ^_ £ i „ 
+ £ ^ Z wL + SZW £ 5Z$l n + 5Z( 2 )] Ao • 

n n / 



(4.5) 



In the equations above, »4j represents the amplitude at i loops stripped of the factor a. In 
the case of the process qq —* ti, the wave function renormalization factors of massless quarks 
vanish at one loop, while the ones of the massive quarks in the on-shell renormalization 
scheme are given by 

«SU = w(£)' '(-s-i^)- (46) 

where the subscript M indicates massive quarks and where C(e) = (47r) e r(l + e). The 
one-loop renormalization constant for a in the MS scheme is given by 

«f'=°Wf^5(-i)- (4 - 7) 

where (3o = 11/6Ca — l/3(iV/ + Nh) and 7 is the Euler-Mascheroni constant 7 0.577216. 
Therefore, the overall one-loop counter term is given by 



52£>„ + 5Z<P= C(£ ' 



WF,M ■ ' ' J . 



2(3 + 3 + 4e + In I At 



m 2 



+ C (e 2 ) . (4. 



To renormalize the two-loop diagrams contributing to the Ni corrections of the par- 
tonic cross section it is necessary to extract from the last two lines in Eq. (|4.5| ) the terms 
proportional to iVj. Taking into account the fact that the wave function renormalization 
factors are zero for the incoming particles and identical for the massive ones, we find that 

10 

4 2 eif i} = A™ + (6zM tM + 2SZ^) A^ + 28Z^J2A^ 

3=3 

+ (szSfj) + 25Z^5Z^ M + 5Z^) Ao . (4.9) 

In Eq. ( [4.9| ), the quantity A^ is the amplitude of the j-th diagram in Fig. ^ (stripped of 
the factor a). The renormalization coefficients not previously defined are: 

6z£' N ^ = C(e) 



r(l +e) 6e ' 
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Figure 4: One-loop diagrams. Thin arrow lines represent massless quarks, thick arrow line massive 
quarks, dashed arrow lines are Faddeev-Popov ghosts, and coiled lines are gluons. 



11 



SZ (2,N t ) = c{e) 2 



e-T* VJVz /5 1 11 



«5Z, 



(2,JVl) 
wp.M 



C( £ ) 



2 / M 



2 \ 2e 



CfNi 



1 



9 59 ?r 2 ^, , 



16e 2 32e 64 24 



(4.10) 



In order to renormalize the part of the squared matrix element proportional to Nh, 
one has to consider the terms proportional to in the last two lines of Eq. ( |4.5| ), and to 
add the counter term for the on-shell mass renormalization: 

10 

A?en h) =4 Nh) + (sZ^Im + 25Z^) A{ d2) + 25Z^ £ A^ ] - 28Z$A^™* CT ] 

3=3 
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+ ( SZSS + + 25Z^5Z^ + «J<VfcA A . (4.11 



It must be observed that in this case also the external massless legs acquire a non van- 
ishing two-loop wave function renormalization factor indicated by SZ^jc^n- The quantity 
^(<2 2 ,mass ct ) i nc ii ca t e s the second diagram in Fig. |3| with a mass counter term insertion in 
one of the internal heavy quark lines. The renormalization constant appearing for the first 
time in Eq. ( [4,ll| ) are: 

SZ (l,N h ) = C() e " 7£ ^ 

sz (i) _ sz w 

g _ 7£ v 2 1 , g 1 n 



«S&S> = OW> (^) " C F N h (j^-^ + O W) . (4.12) 



The renormalization coefficients in Eqs. ( 4.10JL12 ) can be found in [20,46]. 



The renormalization of the functions F( (i = l,h,lh) in Eq. ( |2.5[) is trivial since the 
relevant two-loop diagrams are reducible and involve the insertion of two one-loop fermionic 
vacuum polarization insertions on the gluon propagator in the diagram of Fig. |l[ 



5. Results 



The main result of the present paper is an analytic, non-approximated expression for the 
coefficients Ei,Eh,Di,Dh,Fi, Fih, Fh in Eq. Since such a result is too long to be 

explicitly printed here, we included in the arXiv submission of this work a text file with 
the complete result, which is written in terms of one-dimensional HPLs of maximum weight 
four and two-dimensional HPLs of maximum weight three. Since the coefficients in Eq. ( |2.5| ) 
still contain infrared poles, the result is dependent on the choice of a global, e-dependent 
normalization factor. With our choice, we factor out an overall coefficient 

C 2 (e) = [(^fT(l+e)] 2 . (5.1) 

We also provide two codes, one written in Fortran, the other as a Mathematica package, 
that numerically evaluate the analytic expression of the quantities listed above for arbitrary 
values of the mass scales involved in the calculation. 

In order to cross check our results, we expanded them in the s, \t\, \u\ m 2 limit. The 
first term in the expansion agrees with the results published in [20] ; the second order term 
agrees with the results found in the Mathematica files included in the arXiv version of [24]. 
We also find complete agreement with the numerical result of Table 3 in [24], corresponding 
to a phase space point in which the s, \ t\, \u\ m 2 approximation cannot be applied. 
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It is straightforward to expand our result for values of the center of mass energy close 
to the production threshold. We define 



Am? 



cos 9 



Lfj, = In 



ln 2 = ln(2) 



(5.2) 



where 6 is the scattering angle in the partonic center of mass frame, and we expand our 
results in powers of the heavy quark velocity 0, up to terms of order (5 2 . We find 
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(5.3) 
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F h ((3,0 = ^r + O(P 2 )+O(e), 



(5.8) 




- 4C(2)/3 + O (P 2 ) + O (e) 



(5.9) 



These expansions could be used in the future in the calculation of logarithmically enhanced 
terms at production threshold. 

6. Conclusions and Outlook 

In this paper, we presented the analytic calculation of the two-loop fermionic corrections to 
the heavy-quark production amplitude for qq — ► tt, retaining the exact heavy-quark mass 
dependence. Our work serves as an independent confirmation of recent results obtained 
analytically as small- mass expansions [20] and numerically [24] . We also provide new results 
on the threshold expansion of the amplitude. 

Our result represents a gauge invariant sub-set of the full two-loop virtual correction to 
the partonic process qq — ► tt. In order to complete the analytic calculation of the two-loop 
virtual corrections, it is necessary to calculate the diagrams that do not contain fermion 
loops, which is currently in progress. Likewise, analytic results for the two-loop amplitude 
for gg — > ti could be obtained in the same calculational framework. 

In order to obtain NNLO predictions for the total tt production cross section and for 
differential distributions, it is necessary to combine the two- loop virtual corrections with the 
already available [10] one- loop corrections to the ti+(l parton) process and the tree-level 
ti+(2 partons) process. Since these contribute to infrared-divergent configurations where 
up to two partons can become unresolved, their implementation requires the application of 
a NNLO subtraction method. The methods presently available [47-49] have been applied 
up to now [50-52] to at most 1 — > 3 processes in e + e~ annihilation and 2 — > 1 processes 
at hadron colliders, such that a calculation of a hadronic 2^2 process, involving massive 
partons, will be a new step in complexity, potentially requiring further refinements of these 
methods. 
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A. Harmonic Poly logarithms 

The results of this work are conveniently expressed, in the non-physical region s < 0, 
in terms of one- and two-dimensional HPLs. Nowadays, harmonic poly logarithms are 
extensively used in multiloop computations, therefore in this section we just summarize 
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their definition; the reader interested in the algebraic properties of these functions can find 
detailed discussions of the topic in the available literature [44,45]. 

In the process under study, five different weight functions are needed; they are 



fw(x) 



1 



x — w 



with w € ^ 0,1, —1, —y, — 

y 



The weight-one HPLs are defined as 

G(0;x)=]nx, G(w;x) = [ dtf w (t) 

Jo 

HPLs of higher weight are defined by iterated integrations 



G(w,---;x) = f X dtf w (t)G(- ■ ■ -t) , 
Jo 



(A.l) 



(A.2) 



(A.3) 



with the only exception of the HPLs in which all the weights are zero which are defined as 
follows 



1 



G(0,0,--- ,0;x) = — ln n x 
v ' n! 



(A.4) 



The reader should be aware of the fact that in the original definition of Remiddi and 
Vermaseren, the weight function corresponding to the weight 1 was f± = 1/(1 — x). In 
order to translate the HPLs defined with the Remiddi- Vermaseren convention to the ones 
employed in this work (and vice versa) it is sufficient to multiply each HPL by a factor 
(— l) n , where n is the number of weights equal to 1. 

The weights — y and —1/y were already introduced in [28,45]. In our results, the two- 
dimensional harmonic polylogarithms have maximum weight three. As it is well known, if 
the weight is not larger than three, these functions can be rewritten in terms of Nielsen poly- 
logarithms. For completeness, we list below the explicit expression of the two-dimensional 
harmonic polylogarithms which appear in our analytic results in terms of Nielsen's poly- 
logarithms: 



G(-y;x 
G{-l/y;x 
G(-y,0;x 
G(-l/y,0;x 
G(-y,l;x 

G(-l/y,l;x) 



In 



x + y 
V 



ln(l + xy) , 

Ynix) In f — — - 
V V 



+ Li 2 i 

V V 

In (a;) \n(xy + 1) + Li 2 (— xy) , 
1 



ln 2 (y + 1) - ln(l - x) ln(y + 1) - ln(y) ln(y + 1) 



+ ln(l - x) \n(x + y)- Li 2 (-y) + Li 2 



1 — x 

y + i 



TT 

~6 



- ln 2 (y + 1) - ln(l - x) ln(y + 1) + ln(l - x) \n(xy + 1) 



+Li 2 (-y) +Li 2 



y-xy 

y + i 
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G(-y,0,0;x) = -^ln 2 (x)ln (l + ^ + ln(x)Li 2 ( --) - Li 3 (-- 



2 \ yj \ yJ \ y 

G(-l/y, 0, 0; x) = - ln 2 (x) In (1 + xy) + ln(x)Li 2 (-xy) - Li 3 (-xy) , 

G(-y, 1, 0; x) = -\ ln 3 (l-x)-ln(x) ln 2 (l-x)-ln(y) ln 2 (l-x) + \ ln(y+l) ln 2 (l-x) 

+ ln(x + y) ln 2 (l — x) — — ln 2 (y + 1) ln(l — x) — ln 2 (x + y) ln(l — x) 
— ln(x) ln(y) ln(l — x) — ln(x) ln(y+l) ln(l — x)+ln(y) ln(y+l) ln(l — x) 
+2 ln(x) ln(x + y) ln(l — x) + ln(y) ln(x + y) ln(l — x) 

-Li 2 (x) ln(l - x) + Li 2 f-^ ln(l - x) + Li 2 (-y) ln(l - x) 

+ i ln(x) ln 2 (y + 1) - ^vr 2 ln(x) - ln(x) ln(y) ln(y + 1) - ln(y)Li 2 (x) 

+ ln(x + y)Li 2 (x)-ln(x)Li 2 (-y) +ln(x)Li 2 (j^) + ln (y) Li 2 (^~) 

-ln(x + y)Li 2 f^-)+ln(y)Li 2 filz^) + ln(x + y)Li 2 fc^A 
\x + y) \ x + y J V x + y y 

-ln(y)Li 2 ( I**-) + ln(x + y)Li 2 f ^ - Li 3 (l - x) - Li 3 ^ ' 



x — 1/ \x— ly \x— 1 

"Lis (-£) - Lia(-y) + Li 3 - Li 3 

+Lia f -Li 3 f Ji-) + Li 3 f -Li 3 f ^ + C (3) , 



(x— l)y y \x + y J \ x + y J \x— 1 

G(-l/y, 1, 0; x) = -± ln 3 (l - x) - ln(x) ln 2 (l - x) - X - ln(y) ln 2 (l - x) 

+ - ln(y + 1) ln 2 (l - x) + ln(xy + 1) ln 2 (l - x) - - ln 2 (y + 1) ln(l - x) 

— ln 2 (xy + 1) ln(l — x) — ln(x) ln(y + 1) ln(l — x) 

+21n(x) ln(xy+l) ln(l — x)+ln(y) ln(xy+l) ln(l — x) — Li 2 (x) ln(l — x) 



. / 1 



x 



-Li 2 (— y) ln(l - x) + Li 2 (— xy) ln(l - x) - Li 2 y — — - j ln(l - x) 
-Li 2 ( ^J^ ) ln ( X " x ) + \** ln ( X " x ) + \ ln ( x ) ln2 (y + !) 



+ ln(xy + l)Li 2 (x) + ln(x)Li 2 (-y) + ln(x)Li 2 



y-^y 
y + i 



- ln( * y + 1)Li2 (x^Tl) + ln(a ^ + 1)Li2 (x^Tl) 

+ ln(xy + l)Li 2 (^±1) - Li 3 (1 - x) - Li 3 (-^) + Li 3 (-1 
-Li 3 (-xy) + Li 3 f + Li 3 f ^±il^ - Li 3 ^ " * 



y + ly \ x — 1 y \ y + 1 
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G(-y,0, l;x) = -ln 3 (l-x) + ln(x) ln 2 (l-x) + - ln(y) ln 2 (l-x) -ln(x+y) ln 2 (l-x) 
3 2 

+ ln 2 (x + y) ln(l — x) + ln(x) ln(y) ln(l — x) — ln(x) ln(x + y) ln(l — x) 
- ln(y) ln(x + y) ln(l-x) + Li 2 (x) ln(l-x) + Li 2 ( (1 ~ ^ ) l n (l - x) 



x + y 

+Li 2 (j—f ) ln ^ - z ) - \^ M 1 - x ) + ln(y)Li 2 (x) 

- ln(x + y)Li 2 (x) - ln(y)Li 2 (^^) + + y)Li 2 (^^) 
+ ln(y)Li 2 ( ^) ~ Hx + y)U 2 ( ^) + ln( y )Li 2 (* + y 

V x + y J V x + y J 

- ln(x + y)Li 2 f + Li 3 (l - x) - Li 3 (-y) + Li 3 ( 11 



x — 1 J ' \x + y 

G(-l/y, 0, 1; x) = - ln 3 (l - x) + ln(x) ln 2 (l - x) + - ln(y) ln 2 (l - x) 



- ln(xy + 1) ln 2 (l - x) + ln 2 (xy + 1) ln(l - x) - ln(x) ]n(xy + 1) x 
x ln(l — x) — ln(y) \n(xy + 1) ln(l - x) + Li 2 ^ ~ m (l ~~ x ^ 



1 \ , , w. ( 1- 



+ ln(xy + l)Li 2 - ln(xy + l)Li 2 



xy + iy v * ' ' '"\xy + l 
- ln(xy + l)Li 2 (t^T^) + Li 3 (l - x) - Li 3 

G(-y, 1,1; x) = --ln(y + l)ln 2 (l-x) + -ln(x + y)ln 2 (l-x) 
+Li 2 ( — 4 J ln(l - x) + Li 3 ( — — ) - Li 3 



x 



y + iy \y + 1 J Vy + i 

G(-l/y, 1, 1; x) = Info + 1) ln2 (! " *) + 5 H X V + x ) h 2 (l " x ) 

We first obtained the squared matrix elements in the non-physical region s < 0. The 
corresponding quantities in the physical region s > could be obtained by analytic con- 
tinuation to the complex value s —> s + iS, where 5 —> + . For s > Am 2 the variable x 
becomes 

x = —x + z<5 , (A. 6) 
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where 

x' = ^" V g " 4m2 , (A.7) 
+ V s — Am 2 

So that < x 1 < 1 for 4m 2 < s < oo. While the HPLs with argument y and z are always 
real, the HPLs of x are complex for s > 0. In particular, the imaginary part of the HPLs 
of argument x for s > Am 2 is defined when the analytic continuation of the logarithm is 
specified: 

G(0; x) -> G(0; -x' + %S) = G(0, a/) + m . (A.8) 

In the quantities and -D^ one also encounters the variable x p defined in Eq. ( |3.2j ). In 
this case, in the physical region one finds that 

x p — ► — + z<5 , (A. 9) 

where 



4 = v ; : v ; ; (a.io) 



, V-s + 4m 2 — 
Vs + 4m 2 + \/s ' 

therefore < x p < 1 when < s < oo. The imaginary part of the HPLs of argument is 
defined according to the same principle that determines the imaginary part of the HPLs of 
argument x (Eq,( [A,8| )), The analytic continuation of the two-dimensional HPLs presented 
in this appendix is determined by the analytic properties of the functions In, Li2 and L13. 

B. Master Integrals 

In this Appendix we collect the Mis for the topologies in Fig. §-(k), Fig. §-(1), Fig. ||-(k), 
and in Fig. ||-(1) that are not yet available in the literature. 

The explicit expression of the Mis depends on the chosen normalization of the inte- 
gration measure. The integration on the loop momenta is normalized as follows 

where C(e) was defined in Eq. (pHp. In Eq. pip \i stands for the 't Hooft mass of 
dimensional regularization. The integration measure in Eq. ( |B.l| ) is chosen in such a way 
that the one-loop massive tadpole becomes 

® dk> PT^ 2 = ~4(1 -e)e ' (B ' 2) 

In calculating the squared matrix element, we multiply our bare results by (/U 2 /m 2 ) e , in 
order to make explicit the dependence on the top scale. We also point out that, since 
the squared matrix element still contains soft and collinear divergencies regulated by e, it 
depends on the normalization of the integration measure. In particular, in order to match 
our results with the ones of [20,24], it is necessary to multiply them by the factor 



e 



" 27£ . , 2 



r(l + e) 2 ,v "'~ ' 3 



C(2)e 2 + ^C(3)e d + O (e 4 ) . (B.3) 
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There is a single Mis belonging to the topology Fig. |2j- (k) 




^ d k l r D d k 2 



P (fo + pi) P (k 2 ) Po (h - k 2 ) P m (h + P3 ) ' 



(B.4) 



where we define 



we then find 



P (A:) = k 2 , P m (k) = k 2 + m 2 ; 



(B.5) 




A i e i + 0{e 2 



i=-2 



(B.6) 



A-2 = — , 
32 ' 

A-i = -J-[5y-2(l + y)G(-l;y)] , 
32y 

A = -J- [19y + %C(2) - 10 (1 + y) G (-1; y) 
32y 

+8 (l + y)G (-1, -1; y) - 4G (0, -1; y) - 4yG (0, -1; y)] 



32yL 



65y + 20yC (2) + 8y( (3) - 2 (1 + y) (19 + 4( (2)) G (-1; y) 

+40 (1 + y) G (-1, -1; y) - 20G (0, -1, y) - 20yG (0, -1; y) - 32G (-1, -1, -1; y) 

-32yG (-1, -1, -1; y) + 16G (-1,0, -1, y) + 16yG (-1,0, -1; y) 

+16G (0, -1, -1; y) + 16yG (0, -1, -1; y) -8G (0, 0, -1; y) -'8yG (0, 0, -1; y)] . (B.7) 

Also the four point topology Fig. |2|-(1) has a single MI: 




® d k 1 T> d k 2 



P (h +pi)P (h +pi+ P2) Po (k 2 ) Po (h - k 2 ) P m (h + p 3 ) 



, (B.8) 



with 




1 



(B.9) 



i=-3 



^4-3 
A-2 



32(y + 1) ' 

32 J + 1} (G(0; z) - 2G(1; z) - 2G(-1; y) + 2 ) 
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3C(2) + G(0, 0; x) - 2G(0, 1; x) - 2G(1, 0; x) + 4G(1, 1; x) 



32(y + 1) 

+G(0; x)(2 - 2G(-1; y)) - 4G(-1; y) + G(l; x)(4G(-l; y) - 4) 
+4G(-l,-l;y)-2G(0,-l;i/)+4| , 



An 



1 



2(3C(2) + 4C(3) - 4) + G(0, 0, 0; x) - 2G(0, 0, 1; x) - 2G(0, 1, 0; x) 



32(y + 1) 

+4G(0, 1, 1; x) - 2G (-1/y, 0, 0, x) + 4G (-1/y, 0, 1, x) + 2G (-1/y, 1,0, x) 

-4G (-1/y, 1, 1, x) + 2G(-y, 1,0; x) - 4G(-y, 1, 1; x) + G(l, 1; x)(8 - 16G(-1; y)) 

+G(0, 0; x)(2 - 2G(-1; y)) + 2(3C(2) - 4)G(-1; y) - 2G (-1/y, 0, x) G(-l; y) 

+4G (-1/y, 1, x) G(-l; y) - 2G(-y, 0; x)G(-l; y) + 4G(-y, 1; x)G(-l; y) 

+G(0, 1; x)(4G(-l; y) - 4) + G(l, 0; x)(8G(-l; y) - 4) 

+G(1; x)(2(7C(2) - 4) + 8G(-1; y)) + 8G(-1, -1; y) 

+G(0; x)(-3C(2) - 4G(-1; y) + 4G(-1, -1; y) - 2G(0, -1; y) + 4) 

-4G(0, -1; y) + G(-y; x)(2G(0, -1; y) - 4G(-1, -1; y)) 

-+G (-1/y, x) (-8C(2) - 4G(-1, -1; y) + 2G(0, -1; y)) 

-8G(-l,-l,-l;y)+4G(-l,0,-l;y)+4G(0,-l,-l;y) - 2G(0, 0, -1; y)] . (B.10) 

We now consider the Mis involved in the calculation of the part of the amplitude 
proportional to N^. Topology Fig. |||-(k) has two Mis: 




® d k 1 T) d k 2 



P (h + Pl ) P m (fc 2 ) P m (h - k 2 ) P m (h + P3 ) ' 



(B.ll) 



with 




Ae + 0(e) 



(B.12) 



A_ 2 



32 

-i- [ 5 y-2(y+l)G(-l;y)] 
32y 



1 



$ ^ ^ . [_10G(-1; y)(y + l) 2 + 4G(-1, -1; y)(y + l) 2 - 4G(0, -1; y)(y + 1) 
-4yG(0, 0, -1; y) + y (19y - 4C(3) + 19)] , 

qo J— T s [20G(-1, -1; y)(y + l) 2 - 8G(-1, -1, -1; y)(y + l) 2 
32y(y + 1) L 

+8G(-1, 0, -1; y)(y + 1) + 8G(0, -1, -1; y)(y + 1) - 8(y - 1)G(1, 0, -1; y)(y + 1) 
+4(y - 1)G(1; y)C(2)(y + 1) + 4(y - 2)(2y + 1)G(0, 0, -1; y) - 8yG(-l, 0, 0, -1; y) 
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+8yG(0, -1,0, -1; y) + 8yG(0, 0, -1, -1; y) - 12yG(0, 0, 0, -1; y) 
+16yG(0, 1, 0, -1; y) - 8yG(0, 1; y)((2) + 4G(0, -1; y) (2C(2)y - 5y - 5) 
-2G(-1; y) (2C(2)y 2 + 19y 2 + 4<C(3)y + 38y - 2C(2) + 19) + | (6C(2) 2 + 325y 



+20yC(3) - 40C(3) + 325) 
As second MI for the topology Fig. ^)-(k), we chose 



(B.13) 




D d k 1 r D d k 2 

P 2 (h + Pl ) P m (k 2 ) P m ih - k 2 ) P m (h + Pa) ' 



(B.14) 



where 




_L^ A ^ + G-( £ 2 ) , 



i=-2 



(B.15) 



^4-2 

A 



16(y + l) 

y-i 



16y(y + 1) 
8y(y + l) 2 



G(-l;y), 
(y + l)G(-l;y)-(y + l)G(-l,-l;y) 



2y 



-G(0,-l;y) + -^(y-C(2) + l) 
1-2/ 

8j/( y y + 1 1 )2 [ - 2 (y + - 1 ; y) + 2G (°> -i; i/) + 2 (y + i)Gf(-i, -i, -i; y) 

-2G(-1, 0, -1; y) - 2G(0, -1, -1; y) - 2(y - 1)G(0, 0, -1; y) 

+2(y - 1)G(1, 0, -1; y) - (y - 1)G(1; y)C(2) + G(-l; y) (C(2)y + 2y - C(2) + 2) 



y-i 



(121n(2)C(2)-4C(2)+K(3)-8C(3)) 



(B.16) 



The box topology Fig. ^-(1), has also two Mis: the first one is 



D d k 1 D d k 2 




PO (h +Pl)P (h +P1+ p 2 ) P m (fc 2 ) Pm (h ~ k 2 ) P m (h + P3 ) 



, (B.17) 



with 




(B.18) 



i=-3 
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A- 3 



1 

32(y + 1) ' 
1 



32(x-l)(y + l) 
1 



2G(-l;y)(x - 1) + 2(x - 1) - (1 + x)G(0;x) 



3C(2)x + 4x + 6(x + 1)G(-1, 0; x) + (-5x - 1)G(0, 0; x) 



32(x-l)(y + l) 

-4(x - 1)G(-1; y) + G(0; x)(2(x + 1)G(-1; y) - 2(x + 1)) + 4(x - 1)G(-1, -1; y) 
-2(x-l)G(0,-l;y) + 3C(2)-4 



An 



1 



36(x + l)G(-l,-l,0;x) + 18(x + 1)G(-1, 0, 0; x) 



32(x-l)(y + l) 

+6(5x + 1)G(0, -1, 0; x) + (-5x - 1)G(0, 0, 0; x) - 2(5x - 3)G(1, 0, 0; x) 
+4(x + 1)G(1, 1, 0; x) + 2(x + 1)G (-l/y, 0, 0, x) - 2(x + 1)G (-l/y, 1,0, x) 
-2(x+l)G(-y, 1, 0; x) - 4(x+l)G(l, 0; x)G(-l; y)+2(x+l)G (-l/y, 0, x) G(-l; y) 
+2(x + l)G(-y, 0; x)G(-l; y) + G(-l, 0; x)(12(x + 1) - 12(x + 1)G(-1; y)) 
+G(0, 0; x)(2(5x + 1)G(-1; y) - 2(5x + 1)) + 8(x - 1)G(-1, -1; y) 
-4(x-l)G(0,-l;y) - 2(x + l)G(-y; x)G(0, -1; y) - 8(x - 1)G(-1, -1, -1; y) 
+4(x - 1)G(-1, 0, -1; y) + 4(x - 1)G(0, -1, -1; y) - 2(3x - 1)G(0, 0, -1; y) 
-18(x+l)G(-l;x)C(2) - 4(x+l)G(l;x)C(2) - 2G(-l;y) (3C(2)x+4x+3C(2) - 4) 
+G(0; x)(15C(2)x - 4x + 4(x + 1)G(-1; y) - 4(x + 1)G(-1, -1; y) 
+2(x+l)G(0, -1; y) + 3C(2) - 4) + G (-l/y, x) (2(x+l)G(0, -1; y)+4(x+l)C(2)) 
+2 (3C(2)x + 12C(3)x + 4x + 3C(2) + 4C(3) - 4) 1 . (B.19) 

As second MI for topology Fig. |3|-(1) we chose 



D d k 1 D d k2 



PO (h + pi) P (fa +P1+ P2) Pi (k 2 ) Pm {fa ~ k 2 ) Pm {fa + Ps) 



, (B.20) 



with 




(B.21) 



A-2 
A-! 

A () 



32(x-l)(x + l)(y + l) 
x 

~32(x-l)(x + l)(y + l) 
x 



32(x-l)(x + l)(y + l) 



G(0;x), 

6G(-1, 0; x) - 3G(0, 0; x) + 2G(0; x)G(-l; y) + 3C(2) 
18G(-1, 0, 0; x) - 36G(-1, -1, 0; x) + 18G(0, -1, 0; x) 
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-3G(0, 0, 0; x)-2G(l, 0, 0; x)+4G(l, 1, 0; x) + 2G (-1/y, 0, 0, x) -2G (-1/y, 1, 0, x) 
-2G(-y, 1, 0; x) - 12G(-1, 0; x)G(-l; y) + 6G(0, 0; x)G(-l; y) 
-4G(1, 0; x)G(-l; y) + 2G (-1/y, 0, x) G(-l; y) + 2G(-y, 0; x)G(-l; y) 
-2G(-y; x)G(0, -1; y) - 2G(0, 0, -1; y) - 18G(-1; x)C(2) - 4G(1; x)C(2) 
-6G(-1; y)C(2) + G (-1/y, x) (2G(0, -1; y) + 4<(2)) 

+G(0; x) (-4G(-1, -1; y) + 2G(0, -1; y) + 9C(2)) + 16C(3) ] . (B.22) 
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